from sklearn.svm import SVR
import numpy
import numpy as np
import pandas as pd
from sklearn.metrics import mean_squared_error
from sklearn.metrics import mean_absolute_error
from sklearn.metrics import r2_score
import matplotlib.pyplot as plt


def mean_absolute_percentage_error(y_true, y_pred):
    y_true, y_pred = np.array(y_true), np.array(y_pred)
    return np.mean(np.abs((y_true - y_pred) / y_true)) * 100


dat = pd.read_csv('US-SRG_HH_TOTAL.csv', dtype="float64")
len_dat = dat.shape[0]
# dat = dat[0:len_dat / 3]
# dat = dat[len_dat / 3 - 1:len_dat / 3 * 2]
# dat = dat[len_dat / 3 * 2 - 1:]
# dat = dat[len_dat / 3 - 1:].

max = max(dat.NEE)
min = min(dat.NEE)
dat = (dat - dat.min()) / (dat.max() - dat.min())

dat.dropna(inplace=True)

X = dat.loc[:, [x for x in dat.columns.tolist() if x != 'NEE']].as_matrix()
y = np.array(dat.NEE)
T = 10

train_size = int(X.shape[0] * 0.75)

clf = SVR(gamma='scale', C=1.0, epsilon=0.2)
clf.fit(X[:train_size], y[:train_size])

y_train = clf.predict(X[:train_size - T + 1])
y_train_true = y[:train_size - T + 1]

y_pred = clf.predict(X[train_size:])
y_pred_true = y[train_size:]


def get_r2_numpy(x, y):
    slope, intercept = np.polyfit(x, y, 1)
    r_squared = 1 - (sum((y - (slope * x + intercept)) ** 2) / ((len(y) - 1) * np.var(y, ddof=1)))
    return r_squared


def inverse(data, max, min):
    for i in range(len(data)):
        data[i] *= (max - min)
        data[i] += min
    return data


y_train = inverse(y_train, max, min)
y_train_true = inverse(y_train_true, max, min)
y_pred = inverse(y_pred, max, min)
y_pred_true = inverse(y_pred_true, max, min)

print("train rmse: ", numpy.sqrt(mean_squared_error(y_train, y_train_true)))
print("train mae: ", mean_absolute_error(y_train, y_train_true))
print("train mape: ", mean_absolute_percentage_error(y_train_true, y_train))
print("train r square: ", get_r2_numpy(y_train, y_train_true))
print("train corr: ", np.corrcoef(y_train, y_train_true))

print("predict rmse: ", numpy.sqrt(mean_squared_error(y_pred, y_pred_true)))
print("predict mae: ", mean_absolute_error(y_pred, y_pred_true))
print("predict mape: ", mean_absolute_percentage_error(y_pred_true, y_pred))
print("predict r square: ", get_r2_numpy(y_pred, y_pred_true))
print("predict corr: ", np.corrcoef(y_pred, y_pred_true))

true_list = numpy.array(numpy.concatenate((y_train_true, y_pred_true)))
pred_list = numpy.array(numpy.concatenate((y_train, y_pred)))

# print("true list: " + str(list(true_list)))
# print("pred list: " + str(list(pred_list)))

fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)
ax.plot(range(len(true_list)), true_list, c="black")
ax.plot(range(len(pred_list)), pred_list, c="red")
x_ticks = ax.set_xticks(range(len(true_list) / 24, len(true_list), len(true_list) / 12))
x_labels = ax.set_xticklabels(range(1, 13, 1), fontsize="small")
ax.set_xlabel("Month")
ax.set_ylabel("NEE(gCO2 m-2 d-1)")
ax.legend()
plt.show()
